ohi logo
OHI Science | Citation policy

Trend vs scores

(figure out which regions to label)

index <- read.csv("data/scores_eez2017.csv") %>%
  select(region_id, region_name, Index)

vaxis <- index$Index[index$region_id==0]

trend <- read.csv(sprintf("data/trends_%s.csv", scenario)) %>%
  select(region_id, index_trend = Index) %>%
  left_join(index, by= "region_id") %>%
  filter(region_id != 0)


mod <- lm(index_trend ~ Index, data=trend)
summary(mod)
## 
## Call:
## lm(formula = index_trend ~ Index, data = trend)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.4419 -0.3434  0.0874  0.4413  1.9802 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -1.396047   0.395124  -3.533 0.000501 ***
## Index        0.017490   0.005852   2.989 0.003124 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.7964 on 218 degrees of freedom
## Multiple R-squared:  0.03936,    Adjusted R-squared:  0.03495 
## F-statistic: 8.932 on 1 and 218 DF,  p-value: 0.003124
  p <- ggplot2::ggplot(trend, aes(x=Index, y=index_trend)) +
    geom_point(aes(text=paste0("rgn = ", region_name)), shape=19, size=2, color="gray", alpha = 0.75) +
    geom_point(aes(text=paste0("rgn = ", region_name)), shape=19, size=2, color="gray", alpha = 0.75) +
    geom_hline(yintercept = 0, color="red", size=.5, linetype=3) +
    geom_vline(xintercept = vaxis, color="red", size=.5, linetype=3) +
    stat_smooth(method=lm, se=FALSE, color="black", size=0.5) +
    theme_bw() + 
    labs(y="Average change in OHI score per year)", x="OHI score (2017)") 
## Warning: Ignoring unknown aesthetics: text

## Warning: Ignoring unknown aesthetics: text
  plotly_fig <- plotly::ggplotly(p)
## We recommend that you use the dev version of ggplot2 with `ggplotly()`
## Install it with: `devtools::install_github('hadley/ggplot2')`
  plotly_fig
### adding labels
  labels1 <- c("Bosnia and Herzegovina", "Vietnam", "Philippines",
              "Eritrea", "Equatorial Guinea", "Nicaragua", "Turkey", "Somalia", "Costa Rica",
          "Glorioso Islands", "American Samoa", "Croatia",
              "Estonia", "New Zealand", "Bahamas")
  
  labels2 <- c("South Georgia and the South Sandwich Islands", "Republique du Congo")
  
setdiff(labels1, trend$region_name)
## character(0)
trend <- trend %>%
  mutate(label1 = ifelse(region_name %in% labels1, as.character(region_name), NA)) %>%
  mutate(label2 = ifelse(region_name %in% labels2, as.character(region_name), NA)) 
  

  p <- ggplot2::ggplot(trend, aes(x=Index, y=index_trend)) +
    geom_point(shape=19, size=1.75, color="black", alpha = 0.5) +
    geom_text(aes(label=label1), hjust=0, nudge_x=0.4, size=2.2) + 
    geom_text(aes(label=label2), vjust=0, nudge_y=0.1, nudge_x=-5, size=2.2) +
   geom_point(data=filter(trend, region_name %in% label1), aes(x=Index, y=index_trend), shape=19, size=1.75, color="darkorange") +
   geom_point(data=filter(trend, region_name %in% label2), aes(x=Index, y=index_trend), shape=19, size=1.75, color="darkorange") +
    geom_hline(yintercept = 0, color="red", size=.2, linetype=3) +
    geom_vline(xintercept = vaxis, color="red", size=.2, linetype=3) +
      stat_smooth(method=lm, se=FALSE, color="black", size=0.3) +
    theme_bw() + 
    labs(y="Average change in Index score per year", x="OHI score, 2017") 

plot(p)   
## Warning: Removed 205 rows containing missing values (geom_text).
## Warning: Removed 218 rows containing missing values (geom_text).

ggsave("figures/trend vs score.png", width=7.5, height=5, units=c("in"), dpi=600)
## Warning: Removed 205 rows containing missing values (geom_text).

## Warning: Removed 218 rows containing missing values (geom_text).
ggsave("figures/trend vs score.tiff", width=6, height=4, units=c("in"), dpi=375)
## Warning: Removed 205 rows containing missing values (geom_text).

## Warning: Removed 218 rows containing missing values (geom_text).

pairwise comparison of goal scores

index_2017 <- read.csv("data/scores_eez2017.csv") %>%
  filter(region_id != 0) %>%
  select(-c(region_name, region_id)) %>%
  select(-c(SPP, HAB, ECO, LIV, FIS, MAR, ICO, LSP, Index))

panel.cor <- function(x, y, digits=2, prefix="", cex.cor, ...){
  usr <- par("usr"); on.exit(par(usr))
  par(usr = c(0, 1, 0, 1))
  r <- abs(cor(x, y, use="na.or.complete"))
  txt <- format(c(r, 0.123456789), digits=digits)[1]
  txt <- paste(prefix, txt, sep="")
  if(missing(cex.cor)) cex.cor <- 0.8/strwidth(txt)
  text(0.5, 0.5, txt, cex = (cex.cor * r + .15)*3)
}

pairs(index_2017,
      lower.panel=panel.smooth, upper.panel=panel.cor, 
      pch=16, cex=.5)

## export as png with 700 x 700 dimensions



#seeing ones were significant
# panel.cor <- function(x, y, digits=2, prefix="", cex.cor) 
# {
#   usr <- par("usr"); on.exit(par(usr)) 
#   par(usr = c(0, 1, 0, 1)) 
#   r <- abs(cor(x, y)) 
#   txt <- format(c(r, 0.123456789), digits=digits)[1] 
#   txt <- paste(prefix, txt, sep="") 
#   if(missing(cex.cor)) cex <- 0.8/strwidth(txt) 
#   
#   test <- cor.test(x,y) 
#   # borrowed from printCoefmat
#   Signif <- symnum(test$p.value, corr = FALSE, na = FALSE, 
#                    cutpoints = c(0, 0.001, 0.01, 0.05, 0.1, 1),
#                    symbols = c("***", "**", "*", ".", " ")) 
#   
#   text(0.5, 0.5, txt, cex = cex * r) 
#   text(.8, .8, Signif, cex=cex, col=2) 
# }
# pairs(index_2016,
#       lower.panel=panel.smooth, upper.panel=panel.cor)

Scores vs. population, HDI, and CHI

index <- read.csv("data/scores_eez2017.csv") %>%
  select(rgn_id = region_id, region_name, Index) %>%
  filter(rgn_id != 0) 

hdi <- read.csv("../../../ohiprep/globalprep/supplementary_information/v2017/HDI_data.csv") %>%
  select(rgn_id, hdi=X2015)

# coastal population
pop <- read.csv("../../../ohiprep/globalprep/mar_prs_population/v2017/output/mar_pop_25mi.csv") %>%
  filter(year == 2017) %>%
  select(rgn_id, coastal_pop = popsum) %>%
  mutate(log_coastal_pop = log(coastal_pop + 1))

## chi data
chi <- read.csv("../../../ohiprep/globalprep/supplementary_information/v2016/CHI_data.csv") %>%
  filter(sp_type == "eez") %>%
  select(rgn_id, chi=mean)


data <- index %>%
  left_join(pop, by="rgn_id") %>%
  left_join(chi, by="rgn_id") %>%
  left_join(hdi, by="rgn_id")

data_noNA <- na.omit(data)

mod1 <- lm(Index ~ chi, data=data_noNA)
summary(mod1)
## 
## Call:
## lm(formula = Index ~ chi, data = data_noNA)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -21.305  -5.066   0.884   6.490  17.844 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   74.910      4.066   18.42  < 2e-16 ***
## chi           -2.773      1.012   -2.74  0.00693 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.185 on 143 degrees of freedom
## Multiple R-squared:  0.04988,    Adjusted R-squared:  0.04324 
## F-statistic: 7.507 on 1 and 143 DF,  p-value: 0.006928
mod2 <- lm(Index ~ hdi, data=data_noNA)
summary(mod2)
## 
## Call:
## lm(formula = Index ~ hdi, data = data_noNA)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -20.3787  -4.8493   0.4664   4.9312  14.9979 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   42.048      2.948  14.264  < 2e-16 ***
## hdi           30.559      4.035   7.574 4.12e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 7.094 on 143 degrees of freedom
## Multiple R-squared:  0.2863, Adjusted R-squared:  0.2813 
## F-statistic: 57.36 on 1 and 143 DF,  p-value: 4.12e-12
mod3 <- lm(Index ~ log_coastal_pop, data=data_noNA)
summary(mod3)
## 
## Call:
## lm(formula = Index ~ log_coastal_pop, data = data_noNA)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -22.2407  -4.8944   0.2453   6.0373  17.5409 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      57.4933     5.7716   9.961   <2e-16 ***
## log_coastal_pop   0.4276     0.3810   1.122    0.264    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.36 on 143 degrees of freedom
## Multiple R-squared:  0.008733,   Adjusted R-squared:  0.001801 
## F-statistic:  1.26 on 1 and 143 DF,  p-value: 0.2636
mod4 <- lm(Index ~ chi + hdi, data=data_noNA)
summary(mod4)
## 
## Call:
## lm(formula = Index ~ chi + hdi, data = data_noNA)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -19.4602  -4.8986   0.6099   5.3907  15.7798 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  53.4601     4.3369  12.327  < 2e-16 ***
## chi          -2.9429     0.8453  -3.481 0.000663 ***
## hdi          30.9008     3.8879   7.948 5.33e-13 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.833 on 142 degrees of freedom
## Multiple R-squared:  0.3424, Adjusted R-squared:  0.3331 
## F-statistic: 36.97 on 2 and 142 DF,  p-value: 1.188e-13
mod5 <- lm(Index ~ log_coastal_pop + hdi, data=data_noNA)
summary(mod5)
## 
## Call:
## lm(formula = Index ~ log_coastal_pop + hdi, data = data_noNA)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -20.3655  -4.4935   0.4726   4.9627  15.2772 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      40.0341     5.4418   7.357 1.38e-11 ***
## log_coastal_pop   0.1439     0.3264   0.441     0.66    
## hdi              30.3499     4.0742   7.449 8.35e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 7.114 on 142 degrees of freedom
## Multiple R-squared:  0.2873, Adjusted R-squared:  0.2772 
## F-statistic: 28.62 on 2 and 142 DF,  p-value: 3.615e-11
mod6 <- lm(Index ~ log_coastal_pop + chi, data=data_noNA)
summary(mod6)
## 
## Call:
## lm(formula = Index ~ log_coastal_pop + chi, data = data_noNA)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -21.6531  -4.9675   0.0188   6.0000  19.1289 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      68.7809     7.0169   9.802  < 2e-16 ***
## log_coastal_pop   0.3996     0.3730   1.071  0.28579    
## chi              -2.7434     1.0121  -2.711  0.00754 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.18 on 142 degrees of freedom
## Multiple R-squared:  0.0575, Adjusted R-squared:  0.04422 
## F-statistic: 4.332 on 2 and 142 DF,  p-value: 0.01493
mod7 <- lm(Index ~ log_coastal_pop + chi + hdi, data=data_noNA)
summary(mod7)
## 
## Call:
## lm(formula = Index ~ log_coastal_pop + chi + hdi, data = data_noNA)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -19.4529  -4.8672   0.3955   5.3944  16.1453 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      51.8808     6.2633   8.283 8.37e-14 ***
## log_coastal_pop   0.1103     0.3147   0.350  0.72649    
## chi              -2.9337     0.8484  -3.458  0.00072 ***
## hdi              30.7391     3.9272   7.827 1.08e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.854 on 141 degrees of freedom
## Multiple R-squared:  0.343,  Adjusted R-squared:  0.329 
## F-statistic: 24.54 on 3 and 141 DF,  p-value: 7.786e-13
AIC(mod1, mod2, mod3, mod4, mod5, mod6, mod7)
##      df       AIC
## mod1  3 1025.1346
## mod2  3  983.6492
## mod3  3 1031.2820
## mod4  4  973.7734
## mod5  4  985.4509
## mod6  4 1025.9670
## mod7  5  975.6471
mod8 <- lm(chi ~ log_coastal_pop, data=data)
summary(mod8)
## 
## Call:
## lm(formula = chi ~ log_coastal_pop, data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.37279 -0.34738  0.05498  0.45766  2.73512 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     3.165134   0.131641  24.044  < 2e-16 ***
## log_coastal_pop 0.055023   0.009737   5.651 4.94e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.695 on 219 degrees of freedom
## Multiple R-squared:  0.1273, Adjusted R-squared:  0.1233 
## F-statistic: 31.93 on 1 and 219 DF,  p-value: 4.935e-08
mod9 <- lm(hdi ~ log_coastal_pop, data=data)
summary(mod9)
## 
## Call:
## lm(formula = hdi ~ log_coastal_pop, data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.30878 -0.10859  0.02746  0.11776  0.22842 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     0.575265   0.100805   5.707  6.4e-08 ***
## log_coastal_pop 0.009349   0.006654   1.405    0.162    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.146 on 143 degrees of freedom
##   (76 observations deleted due to missingness)
## Multiple R-squared:  0.01361,    Adjusted R-squared:  0.006716 
## F-statistic: 1.974 on 1 and 143 DF,  p-value: 0.1622
p <- ggplot2::ggplot(data, aes(y=Index, x=hdi)) +
    geom_point(shape=19, size=2, color="gray", alpha = 0.5) +
    theme_bw()  + 
    stat_smooth(method=lm, se=FALSE, color="black", size = 0.5)+
    labs(x="Human Development Index", y="OHI index score") 

plot(p)   
## Warning: Removed 76 rows containing non-finite values (stat_smooth).
## Warning: Removed 76 rows containing missing values (geom_point).

ggsave("figures/OHIvsHDI.png", height=4, width=4.5)
## Warning: Removed 76 rows containing non-finite values (stat_smooth).

## Warning: Removed 76 rows containing missing values (geom_point).
p <- ggplot2::ggplot(data, aes(y=Index, x=chi)) +
    geom_point(shape=19, size=2, color="gray", alpha = 0.5) +
    theme_bw()  + 
    stat_smooth(method=lm, se=FALSE, color="black", size = 0.5)+
    labs(x="Cumulative Human Impact", y="OHI index score") 

plot(p)   

ggsave("figures/OHIvsCHI.png", height=4, width=4.5)

p <- ggplot2::ggplot(data, aes(y=Index, x=log_coastal_pop)) +
    geom_point(shape=19, size=2, color="gray", alpha = 0.5) +
    theme_bw()  + 
    stat_smooth(method=lm, se=FALSE, color="black", size = 0.5)+
    labs(x="Coastal population (ln)", y="OHI index score") 

plot(p)   

ggsave("figures/OHIvsCoastalPop.png", height=4, width=4.5)

p <- ggplot2::ggplot(data, aes(x=chi, y=log_coastal_pop)) +
    geom_point(shape=19, size=2, color="gray", alpha = 0.5) +
    theme_bw()  + 
    stat_smooth(method=lm, se=FALSE, color="black", size = 0.5)+
    labs(x="Coastal population (ln)", y="Cumulative Human Impact") 

plot(p)   

ggsave("figures/CHIvsCoastalPop.png", height=4, width=4.5)

Trend analysis

### prepare data:
data <- read.csv(sprintf('../OHI_final_formatted_scores_%s.csv', saveDate)) 

data <- data %>%
  filter(dimension == "score") %>%   # focus only on score data
  filter(region_id != 0) %>%         # this weighted mean includes high seas and Antarctica
  filter(region_id <= 250) %>%       # get rid of high seas regions
  filter(region_id != 213) %>%
  filter(!is.na(value))
data_lm <- data %>%
  group_by(goal, region_id) %>%
  do(mdl = lm(value ~ scenario, data = .))

#data.frame(glance(data_lm, mdl))
results_lm <- tidy(data_lm, mdl) %>%
  ungroup() %>%
  filter(term == "scenario") %>%
  arrange(goal, estimate) %>%
  left_join(goal_names) %>%
  left_join(rgn_names) %>%
  select(goal=goal, goal_long=long_goal, region_name = country, region_id, average_change_per_year = estimate, p.value) %>%
  data.frame()
## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable

## Warning in summary.lm(x): essentially perfect fit: summary may be
## unreliable
## Joining, by = "goal"
## Joining, by = "region_id"
# 
# weights <- read.csv("../../eez2016/layers/rgn_area.csv") %>%
#   select(region_id=rgn_id, area_km2)
# 
# results_lm <- results_lm %>%
#   left_join(weights, by="region_id")
# 
# summary_lm <- results_lm %>%
#   group_by(goal, goal_long) %>%
#   summarize(mean=mean(average_change_per_year),
#             mean_wt=weighted.mean(average_change_per_year, w=area_km2))


# color coded histogram of trends

trends_hist <- results_lm %>%
  filter(goal=="Index") %>%
  arrange(average_change_per_year)

values <-   brewer.pal(11, "RdYlBu")
values <- colorRampPalette(values, space = 'Lab')(13)[2:13]
col.brks  <- c(-100, -5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 100)

trends_hist$colors <- cut(trends_hist$average_change_per_year, breaks=col.brks, include.lowest=TRUE, labels=values)
cols <- as.character(unique(trends_hist$colors))


g <- ggplot(trends_hist, aes(x=average_change_per_year)) +
geom_histogram(aes(fill=colors), bins=25, center=1.01, col="gray") +
scale_fill_manual("", breaks=colors, values=cols, guide=FALSE) +
labs(y="Number of regions", x="Average yearly change") +
#geom_vline(xintercept=0, color="red") +
theme_bw()

print(g)

 ggsave("figures/trend_histogram_color.tiff", width=4, height=3, units=c("in"), dpi=300)
data_scores <- data %>%
  filter(goal=="Index") %>%
  filter(scenario == 2017) 

reds <-  colorRampPalette(c("#A50026", "#D73027", "#F46D43", "#FDAE61", "#FEE090"), space="Lab")(66)
blues <-  colorRampPalette(c("#E0F3F8", "#ABD9E9", "#74ADD1", "#4575B4", "#313695"), space="Lab")(35)
colors <-   c(reds, blues)

col.brks <- seq(0,100, by=1)
col.assign <- cut(col.brks, breaks=col.brks, include.lowest = TRUE)
colors_df <- data.frame(colors_plot=colors, col.assign)
colors_df$col.assign <- as.character(colors_df$col.assign) 

data_scores$col.assign <- cut(data_scores$value, breaks=col.brks, include.lowest=TRUE) 
data_scores$col.assign <- as.character(data_scores$col.assign)

data_scores <- left_join(data_scores, colors_df, by="col.assign") %>%
  mutate(colors_plot = as.character(colors_plot)) %>%
  arrange(value)

cols <- as.character(unique(data_scores$colors_plot))

g <- ggplot(data_scores, aes(x=value)) +
#geom_histogram(aes(fill=colors), bins=25, center=1.01, col="gray") +
geom_histogram(aes(fill=col.assign), bins=31) +  
scale_fill_manual("", values=cols, guide=FALSE) +
xlim(c(0, 100)) +
labs(y="Number of regions", x="Index score") +
theme_bw()

print(g)

 ggsave("figures/score_histogram_color.tiff", width=4, height=3, units=c("in"), dpi=300)
 ggsave("figures/score_histogram_color.png", width=4, height=3, units=c("in"), dpi=300)